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I.  INTRODUCTION 

When  considering  dynamical  systems  with  isolated  feedback  mechanisms  or  coupling  devices  to  connect  a  network, 
there  always  exists  a  finite  time  for  the  signal  to  be  transmitted  from  one  device  to  another.  Modeling  such  finite 
transmission  times  requires  the  inclusion  of  one  or  more  delays  in  the  dynamics.  Evidence  of  the  pervasiveness  of 
delays  can  be  seen  in  several  reviews,  such  as  the  application  of  delay  differential  equations  in  many  fields  [1] ,  where 
the  effects  of  delay  are  seen  in  physics,  biology,  and  engineering.  The  analysis  of  the  effects  of  delay  on  stability 
and  bifurcations  is  seen  in  several  collections,  such  as  [2-4].  Other  areas  involving  delay  occur  in  models  involving 
populations  such  as  epidemics  [5-7],  engineering  such  is  machine  tool  cutting  [8],  and  control  [9]. 

Major  applications  and  the  effects  of  delay  and  noise  on  dynamical  systems  can  seen  in  experiments  and  models 
from  nonlinear  optics,  such  as  laser  self- feedback  [10],  delay  coupled  semi-conductor  and  solid  state  lasers  [11-13], 
and  fiber  lasers  [14].  In  many  of  the  coupled  cases,  zero  -lag,  also  called  isochronal,  synchronization  is  an  unstable 
state  [12].  Combining  delayed  self-feedback  along  with  delay  coupling  stabilizes  the  isochronal  state  in  lasers  [15], 
even  in  the  presence  of  noise. 

From  a  modeling  perspective,  the  interaction  of  differential  delay  with  noise  leads  one  to  study  stochastic  effects 
on  infinite  dimensional  systems.  In  approximating  the  phase  space,  finite  dimensional  formulations  based  on  coupled 
map  lattices  were  use  to  initially  examine  approximations  to  steady  state  densities  [16].  Here  the  authors  employ 
a  high  dimensional  Perron-Frobenius  formulation  to  extract  the  density  in  the  thermodynamic  limit.  However,  the 
drawbacks  to  this  formulation  come  from  the  high  dimensionality  of  the  problem  making  analysis  of  the  distribution 
exponents  difficult,  and  the  restriction  to  the  type  of  noise  that  may  be  considered.  For  small  delays,  the  stochastic 
delay  problem  for  Wiener  processes  has  been  re-formulated  using  a  Fokker-Planck  (FP)  approach  [17],  where  the 
new  FP  equations  contain  a  conditionally  averaged  drift  term  in  the  approximation.  The  multivariate  version  of  the 
perturbed  problem  is  derived  in  general  settings  in  [18,  19]  to  both  Ito  and  Stratonovich  formulations. 

Noise  induced  switching  was  explored  numerically  in  bistable  potentials  and  has  been  considered  in  the  presence  of 
delayed  feedback  for  systems  modeling  lasers  to  illustrate  resonant  behavior  due  to  delayed  feedback  [20] .  Stochastic 
resonance  in  bistable  systems  with  delay  and  correlated  noise  sources  was  also  studied  using  a  perturbed  Fokker- 
Planck  approach  for  small  delay [21].  In  addition,  others  have  investigated  noise  induced  switching  in  bistable  systems 
with  delay,  [21-23].  Overall,  many  of  the  above-mentioned  techniques  use  a  perturbed  FP  approach  to  describe  small 
delay  effects  on  the  probability  density  function. 

Rates  of  escape  predicted  from  switching  between  co-existing  stable  states  in  the  absence  of  delayed  feedback  is 
well  investigated  for  a  large  variety  of  well-controlled  systems  and  oscillators  [24-31].  Fluctuations  in  these  systems 
are  usually  due  to  thermal  or  externally  applied  Gaussian  noise.  An  increasingly  important  role  may  be  played  also 
by  more  general  noise  sources,  which  typically  are  non-Gaussian.  It  may  arise,  for  example,  from  just  a  few  two-state 
discrete  particles  jumping  at  random  between  the  states,  often  characterized  as  telegraph  noise  [32].  To  reveal  how 
switching  in  the  presence  of  delay  with  general  noise  sources,  we  consider  a  different  formulation  from  the  usual 
perturbation  procedure  of  the  conditional  averaged  drift  terms. 

Here  we  consider  predicting  the  escape  rates  from  of  a  basin  of  attraction  for  stochastic  systems  with  delayed 
feedback.  The  approach  we  take  is  to  formulate  the  path  which  optimizes  the  probability  of  escape  from  a  stable 
state  with  a  basin  boundary.  The  path  which  optimizes  the  probability  of  escape,  the  optimal  path,  will  be  found 
for  general  noise  sources  using  a  variational  formulation,  originally  posed  for  driven  systems  in  [33,  34],  and  recently 
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used  to  discuss  Poisson  induced  switching  in  [35,  36]. 

The  layout  of  the  paper  is  such  that  in  Sec.  II  we  use  a  variational  approach  to  determine  the  switching  rate  out 
of  the  basin  of  the  attracting  steady  state.  In  Sec.  Ill  we  apply  our  theory  to  a  specific  example  with  additive  noise, 
while  in  Sec.  IV  we  apply  the  theory  to  the  case  of  multiplicative  noise.  In  both  cases  we  obtain  excellent  fit  between 
numerical  simulations  and  the  theory.  We  finish  the  paper  with  a  discussion  in  Sec.  V. 


II.  GENERAL  FORMULATION 

We  consider  switching  processes  for  the  n-dimensional  delay  system  posed  as  a  Langevin  problem: 

x(t)  =  F(x(t),  x(t  -  t))  +  g(x(t))£(t),  (1) 

where  £(£)  is  a  noise  vector  with  independent  components,  and  t/(x(£))  =  diag{g\{t),  g2(t), . . . ,  gn(t)},  a  diagonal 
matrix,  and  r  >  0.  Since  we  are  studying  switching  from  a  given  state,  we  assume  that  in  the  absence  of  noise  that 
there  exists  an  attracting  steady  state,  xa,  and  a  saddle  point,  xs,  where  the  latter  is  on  the  basin  boundary  of  xa; 
i.e.,  F(xa,  xa)  =  F(xs,  xs)  =  0. 

General  noise  characterization  of  £  may  be  accomplished  through  the  characteristic  functional  [37],  where 

(2) 

and  where  (.  .  .}$  means  averaging  over  £(£).  For  Gaussian  noise,  if  D  denotes  the  noise  intensity,  we  characterize 
£(£)  by  its  probability  density  functional  ^[^(t)]  =  exp  (—7 Z^/D), 

=  \f  dtdt'm^t-t'W),  (3) 

where  J~(l  —  t')/2D  is  the  inverse  of  the  pair  correlator  of  £(i).  We  also  assume  the  noise  intensity  D  is  sufficiently 
small  so  that  in  our  analysis  sample  paths  will  limit  on  an  optimal  path  as  D  — >  0. 

If  we  examine  the  tail  of  the  distribution  for  a  large  fluctuation,  assumed  to  be  a  rare  event,  then  the  probability 
of  a  large  fluctuation  is 

■Px[x]  =  exp (-R/D),  R  =  min7?.[x,  £,  A],  (4) 

where 

77[x,£,A]  =  77s  [£(£)]  +  J  A(f)[x(f)  -  F(x,  xT)  -  g(x(t))£(f)]di,  (5) 

and  xr  =  x(f  —  r).  The  approach  we  take  here  is  one  of  a  Hamiltonian  method,  with  added  constraints  [33,  34].  An 
alternative  method  is  to  use  a  Lagrangian  description,  but  then  one  is  forced  to  invoke  an  inverse  function  theorem 
to  convert  it  a  Hamiltonian.  In  many  cases,  as  in  general  descriptions  of  finite  population  interactions  [38],  it  is 
not  possible  to  use  the  Lagrangian  approach.  In  the  general  formulation  above,  the  Lagrange  multipliers,  A(f), 
correspond  to  the  conjugate  momenta  in  the  Hamiltonian  derivation.  They  will  turn  out  to  be  a  deterministic  model 
of  the  effective  force  of  noise  on  the  system,  which  drives  dynamics  towards  the  saddle,  xs,  from  the  attracting  state, 
xa-  Once  sufficiently  past  the  saddle,  the  deterministic  dynamics  takes  over  and  the  process  of  escape  is  completed. 

A.  Equations  of  motion  describing  the  optimal  path 

From  here  on,  we  concentrate  on  the  case  where  we  have  a  scalar  state  variable,  x,  and  the  noise  is  delta  correlated 
in  time.  However,  it  is  clear  that  the  resulting  equations  of  motion  derived  may  easily  be  extended  to  the  general 


Vi  [k]  =  (  exp 


i  J  dfk(t)£(f) 
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vector  case.  To  determine  the  exponent,  7Z,  we  seek  the  equations  which  describe  the  maximum  probability  of  the 
reaching  the  state  at  xs  if  the  initial  state  is  xa,  where  the  exponent  to  extremize  is  now  given  by 


H[x,  £,  A]  =  ^  J  £2{t)dt  +  J  A (t)[x(t)  -  F(x,xT)  -  g(x(t))£(t)\dt. 


(6) 


We  derive  the  variation  (<577.)  by  varying  deviations  from  the  path  that  minimizes  77.  We  start  with  the  variation 
with  respect  to  the  noise  (£(i)). 

S1Z 


=  K[x,t  +  ri,  A]  -77[z,£,A] 


\  j  [(£  +  v)2  -  C2]  dt  -  j  A  [g(Z  +  ri)  -  g£]  dt 
/« 


=  J  (£-  Ag)?7  dt  +  0(r/2). 

Since  g  is  an  arbitrary  smooth  function,  we  have 

£  =  Ag. 

Using  similar  arguments,  we  find  the  variation  with  respect  to  the  Lagrange  multiplier  (A(f)): 

x  =  F  +  g£. 


(7) 

(8) 
(9) 

(10) 

(11) 


Finally,  we  find  the  variation  with  respect  to  x.  Here  the  delay  is  crucial  in  that  it  will  introduce  an  advanced 
term  in  the  equations  of  motion.  The  variation  is  given  as: 

^  =  7Z[x  +  77, 5,  A]  —  72.  [x,  £,  A] 

=  J  A(f)  [77  -  {F(x  +  77,  xT  +  v)  -  F(x,  xT)  -  [g{x  +  rj)  -  g(a:)]£)]  dt 

J 


=  /  A(t) 


v- 


SF(x, 


-t)  - 


5F(x,  xT 


v- 


Sg(x,x  r 


dt  +  Oirj2). 


5x  '  Sx-r  '  8x 
Using  integration  by  parts  in  the  usual  way,  and  rescaling  the  limits  we  derive  the  equation  of  motion: 

-A  -  X(t)^-(x,xT)  -  \(t  +  t)^-{x-t,x)  -  A (t)^(x)£(t)  =  0. 

The  three  coupled  equations  derived  above  may  be  rewritten  as  a  second  order  system: 

x  =  A  g2(x)  +  F(x,Xt) 

A  =  -\2g(x)^(x)  - \^(x,xT)  -  A-t^-(x-t,x), 


which  has  Hamiltonian 


H(: r,  iT,  A)  = 


A 2g2(x) 


+  A  F(x,xT). 


(12) 

(13) 

(14) 

(15) 

(16) 


Notice  the  appearance  of  advanced  terms  in  the  equations  of  motion.  The  equations  of  motion  in  this  formulation 
are  given  by: 


dH  (  \ 

Xo  =  -q^{X0,P0,X0t) 

dH  .  ,  dH  ,  .  .  .  .  /  \  \ 

Po  =  --^{Xo,PQ,X0t)  -  Q^(Xo(t  +  T),po(t  +  T),Xo{t)). 


(17) 


An  alternative  derivation  of  Eq.  17  is  presented  in  the  appendix  for  arbitrary  dimensions.  When  using  the  small 
noise  techniques  for  n-dimensional  stochastic  systems  without  delay,  the  additional  conjugate  momenta  doubles  the 
phase  space  to  2 n.  Attractors  become  saddles,  and  saddles  remain  saddles  but  in  2n-dimensional  space.  If  delay  is 
now  included,  we  discover  the  equations  of  motion  contain  both  retarded  and  advanced  terms  both  of  magnitude  r, 
effectively  doubling  the  time  history  to  2 r  as  well  as  the  phase  space  dimension  to  2 n. 
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B.  Boundary  Conditions 


Recall  that  in  the  absence  of  noise  (£(t)  =  0),  F  has  two  steady  states:  an  attractor  and  a  saddle;  i.e.,  F(xa,  xa)  = 
F(xs,xs)  =  0.  Taking  the  limit  as  t  — >  oo,  we  assume  x(t)  — >  xa,  £(t)  — >  0,  and  A (t)  — >  0.  For  escape  from  the  basin 
of  xa,  as  t  — »  oo,  we  assume  x(t)  — >  xs,  (, {t )  — >  0,  and  A (t)  — >  0. 

To  verify  the  conjugate  momenta  have  the  right  limits,  we  linearize  Eq.  (15)  about  steady  state  values: 

-A  -  A (t)^-(x,x)  -  A (t  +  t)  ^-(x,x)  =  0,  (18) 

where  x  is  either  xa  or  xs ■  Assuming  A  =  eat, 

dF  dF 

-a  - —(x,x)  -  ea'r -^-(x,x)  =  0.  (19) 

We  impose  sufficient  conditions  on  F  so  that  the  linearized  solution  has  the  correct  asymptotic  limits  when  the  delay 
is  small.  More  specifically,  for  the  examples  we  consider  in  the  following  sections,  we  use  Eq.  (19)  to  compute  the 
growth  constant  a  and  confirm  that  we  obtain  the  correct  limiting  behavior  under  the  influence  of  noise. 


C.  Computing  the  exponent  along  the  optimal  path 


With  the  assumptions  above  for  the  equations  of  motion  and  asymptotic  boundary  conditions,  computing  the 
action  7 Z  requires  the  solution  to  the  two  point  boundary  value  problem.  We  suppose  that  the  solution  exists  for 
the  problem  when  r  =  0,  and  that  solutions  for  non-zero  delay  remain  close  to  the  zero  delay  solution.  We  do  not 
necessarily  need  to  assume  the  delay  is  small,  as  done  in  [19].  Then  we  may  compute  perturbations  to  the  action. 

In  order  to  compute  the  correction  to  the  action,  assume  that  perturbations  to  the  optimal  path  when  r  ^  0 
remain  small.  That  is,  if  x(t)  is  the  solution  when  r  =  0,  then  the  perturbation  STx(t)  =  x(t)  —  x(t  —  r)  should 
remain  small.  We  will  need  to  put  conditions  on  the  growth  of  a  general  nonlinear  function,  F,  so  also  assume 
| F{x,  x)  —  F{ x,  xT)\  ;$  k\x  —  xT |. 

The  action  may  now  be  expressed  as  the  following  perturbation  problem: 

ll[x,£,X]  =  1Zo[x,£,  A]  +  1Zi[x,£,  A],  (20) 

where 

Fo[x,  £,  A]  =  ^  J  f(t)dt  +  j  A(t)[£(t)  -  F(x,x)  -  g(x(t))£{t)]dt  (21) 

Ri[i,(,A]  =  J  X(t)[F(x,  x)  —  F(x,  xT)]dt.  (22) 

The  minimizing  solution  is  solved  by  first  considering  the  equations  that  minimize  IZq,  denoted  by  [xa,  C,  A0].  The 
first  order  correction  is  then  evaluated  at  the  zeroth  order  solution.  The  minimizing  solutions  satisfy 

x0  =  X0g2(x0)  +  F(x0,Xo)  (23) 

(24) 


Ao  —  A og(Xo')  q  (*^o)  q  (Xo^Xo*)  X°  ^  (Xo>Xo)- 


The  general  solution  to  Eq.  (21)  will  have  the  form 


A 0(t)  =  —2 ’  and  =  —F(x0(t),Xo(t)). 


(25) 


g2(x0(t)) 

In  the  next  two  sections  we  will  consider  a  specific  example  of  F  for  both  additive  (g(x)  =  1)  and  multiplicative 
{g(x)  =  y/x)  noise.  Using  Eqs.  (23)-(25)  we  compute  the  action,  1Z,  and  then  the  switching  rate  out  of  the  attracting 
basin. 
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III.  ADDITIVE  NOISE  SWITCHING 

A  natural  place  to  discuss  the  escape  rate  is  to  consider  additive  noise;  i.e.,  g(x)  =  1.  For  additive  noise,  either 
the  Hamiltonian  or  Lagrangian  formulation  will  work  to  describe  the  trajectory  of  the  switching  optimal  path.  For 
additive  noise,  the  second  order  delay  system  in  Eq.  (14)  is  equivalent  to  the  equations  of  motion  derived  from  the 
Lagrangian,  L,  which  is  presented  in  the  appendix.  In  the  following  example,  we  use  the  variational  formulation  to 
find  the  escape  rate  exponents  for  a  scalar  delay  example. 


A.  Analysis  for  a  given  function  F: 


As  an  example,  consider  F( x,xT)  =  x(l  —  x)  —  737,  where  the  steady  states  are  xa  =  1  —  7,  xs  =  0,  and 

0  <  7  <  1.  We  first  analyze  the  stability  near  equilibria  to  verify  that  A (t)  approaches  0  as  t  — >  ±oo.  Let  A (t)  be  the 

linearized  solution  to  Eq.  (19)  about  the  steady  states.  Assuming  A(t)  =  exp(at),  for  small  r  we  find  that  near  xa, 
a  =  IZr-y  >  0,  implying  A(t)  — >  0  as  t  — >  — oo.  On  the  other  hand,  in  a  neighborhood  of  xs  we  find  a  =  —  <  0, 

and  A(f)  — >  0  as  t  — >  oo. 

To  compute  the  optimal  path  from  xa  to  xs  that  minimizes  IZo,  we  use  Eq.  (25),  which  for  g(x)  =  1  are 
x0(t)  =  —F(x0(t),  x0(t)),  and  A0(t)  =  —2F(x0(t),x0(t)).  From  Eq.  (10)  we  also  have  that  the  optimal  noise  is 
given  by  £0(t)  =  A 0(t).  For  the  example  F(x,x T)  above,  IZo  has  as  its  solution  x0{t)  =  xa/{  1  +  eXAt).  Notice  that 
limt-Kx,  x0(t )  =  0  =  xs  and  limt^_oo  x0(t)  =  1  —  7  =  xa-  In  addition,  it  easy  to  show  TZo  =  x3A/2>. 

Using  F(x,  x)  —  F(x,  xT)  =  —7(1  —  xT),  we  now  find  the  first  order  action, 


111  =  -27 


/OO 

Xo 

-00 


—  Xo{t  —  r))dt. 


Substituting  for  xo,  the  first  order  results  for  the  action  are  now  given  by, 

^  ..  0.7  ,  x\eXAT{eXAT  -1-xat) 

Ul  ~  ~27  [~~2  + - (1  -  e^r)2 - 


The  full  action  to  first  order  in  5x(t)  is  then 


U(t) 


x3a 


-27  -^  + 


x2a  x2AeXAT  (eXAT -1-xat) 


(1  -  e^r)2 


Recall  xa  =  1  —  7-  Asymptotically,  the  expression  for  small  r  becomes: 

(1-7)3 


72.(7 


-(1  -77)  +  0(t  ). 


(26) 


(27) 


(28) 


(29) 


This  is  our  relationship  for  the  action  on  the  delay. 


B.  Numerical  Comparisons 

We  now  compare  the  theoretical  results  with  numerical  simulations.  In  particular,  we  have  computed  the  mean 
switching  time  as  the  noise  drives  the  system  out  of  the  basin  of  attraction  of  xa  to  x  <  xs-  In  each  of  Figs.  1-3  the 
solid  curve  represents  the  theory,  while  the  data  points  are  the  mean  values  taken  over  1000  simulations.  To  insure 
that  trajectory  doesn’t  drift  back  to  the  original  basin  we  require  that  the  trajectory  be  driven  well  past  xs-  More 
specifically,  we  define  the  switching  time  to  be  when  x  <  —0.2  and  find  that  this  threshold  gives  consistently  good 
results. 
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The  switching  rate  is  proportional  to  the  probability  of  large  fluctuations  "Px[x],  Eq.  (4),  and  is  given  by 

Ws  =  cexp(—R/D),  (30) 

where  for  R  we  use  the  minimum  value  given  by  Eq.  (29).  The  switching  time  is  simply  the  inverse  of  the  rate  and 
so  Ts  =  1/W.  When  plotting  Ts  on  the  log10  scale,  the  constant  c  corresponds  to  a  vertical  shift  in  the  theoretical 
result.  Because  we  are  using  Gaussian  noise,  we  compute  the  prefactor  c  based  on  Kramer’s  theory  [39]  and  find 
that  it  is  given  by  c  =  27r(l  —  7t)/(1  —  7). 

In  Fig.  1,  we  plot  the  switching  time  as  a  function  of  the  inverse  of  the  noise  intensity.  We  see  that,  as  expected, 
when  the  noise  intensity  is  increased  the  switching  time  decreases  and  we  have  excellent  fit  between  the  theory  and 
simulations  over  four  orders  of  magnitude.  In  the  top  curve  the  delay  is  set  to  zero  so  that  F(x,  xT)  =  F(x,  x),  while 
in  the  bottom  curve  the  delay  is  r  =  0.5.  We  see  that  even  though  r  is  not  appreciably  small  we  still  have  excellent 
fit  over  the  entire  range  of  the  noise  intensity.  Further,  the  effect  of  the  delay  is  to  decrease  the  switching  time  or, 
in  other  words,  the  delay  increases  the  probability  of  switching. 

In  Fig.  2  we  fix  the  noise  intensity  at  two  different  values  and  vary  the  delay  r.  As  mentioned  above,  increasing 
the  delay  decreases  the  switching  time  and  our  theory  based  on  small  r  maintains  good  fidelity  up  almost  r  =  1;  for 
r  >  1  the  theory  and  simulation  results  begin  to  diverge. 

In  Fig.  3  the  delay  is  fixed  but  the  strength  of  the  delay  term  7  is  varied  from  zero  to  almost  one.  Recall  that 
xa  =  1  —  7  so  that  increasing  7  moves  the  attractor  closer  saddle  steady  state.  The  effect  is  that  the  noise  has  to 
drive  the  system  a  shorter  distance  before  switching  occurs  and,  hence,  the  switching  time  decreases.  We  note  that 
the  theoretical  curve  deviates  from  the  simulation  results  as  7  gets  close  to  1;  this  is  due  to  a  singularity  in  the 
expression  for  the  shift  constant  c. 

Finally,  the  results  for  the  switching  time  as  a  function  of  7  provide  insight  into  why,  for  this  system,  the  delay 


Figure  1:  Log  of  escape  time  vs.  inverse  noise  intensity  (7  =  0.1).  The  solid  curve  represents  the  theory 
from  Eqs.  (29)-(30),  while  the  data  points  are  the  mean  values  taken  over  1000  simulations. 
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Figure  2:  Log  of  escape  time  vs.  delay  (7  =  0.1).  The  solid  curve  represents  the  theory  from  Eqs.  (29)-(30), 
while  the  data  points  are  the  mean  values  taken  over  1000  simulations. 


Figure  3:  Log  of  escape  time  vs.  dissipation  (r  =  0.1).  The  solid  curve  represents  the  theory  from  Eqs.  (29)- 
(30),  while  the  data  points  are  the  mean  values  taken  over  1000  simulations. 

decreases  the  switching  time.  For  small  x  the  deterministic  repulsive  force  of  the  saddle  is 

x  ~  x  —  7®(t  —  r).  (31) 

For  larger  t  the  second  term  takes  values  of  x  further  back  in  the  history  of  x,  which  because  x  is  generally  decreasing 
from  xa  to  xs  =  0,  are  larger.  The  second  term  being  larger  decreases  the  repulsive  force  of  the  saddle.  Thus,  as  we 
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Figure  4:  Log  of  escape  time  vs.  inverse  noise  intensity  (7  =  0.1).  The  solid  curve  represents  the  theory 
from  Eqs.  (30)  and  (35),  while  the  data  points  are  the  mean  values  taken  over  2000  simulations. 

have  found,  the  delay  makes  it  easier  for  the  noise  to  drive  the  system  out  of  the  basin  of  xa  and  the  switching  time 
decreases. 


IV.  MULTIPLICATIVE  NOISE  SWITCHING 


We  now  consider  the  switching  rate  for  multiplicative  noise  to  illustrate  the  generality  of  the  variational  ap¬ 
proach.  Let  the  multiplicative  function  be  g(x)  =  yT.  The  path  that  minimizes  Eq.  (21)  satisfies  A 0(t)  = 
— 2F(x0{t),  x0(t))/g2(x0(t))  and  x0(t)  =  —F(x0(i),x0(t)).  For  the  example  F(x,  xT)  above,  TZo  has  as  its  solution 
x0(t)  =  xa/{X  +  eXAt)  and  A0(t)  =  —2(1  —  xo  —  7).  We  evaluate 

2  1  roo  r  orn/__  „\l2 


TZo  = 


E( 


2  ,/_^o  V  \/ Xo 

Using  F( x,  x)  —  F(x,x-r)  =  —7(1  —  xT),  we  now  find  the  first  order  action, 

x0(t)[x0(t)  —  x0  (t  —  t)] 


2 F(xo,x0)\  dt  =  x2A 

Xo 


IZi  =  -27 


f 


Xo{t) 


dt 


After  substituting  for  xo  and  evaluating  the  integrals  for  7 Zi  and  TZoo,  we  find  that 


7 Z(r)  «  xa  -  27  (  xa  + 


x2ATeXAT 
—  1  4-  eXAT 

Asymptotically,  the  expression  for  small  r  becomes: 

7 Z(t)  fts  x\{l  -  7 r)  +  0(t2). 


(32) 


(33) 


(34) 


(35) 


Again,  we  compare  the  theoretical  results  with  numerical  simulations.  For  multiplicative  noise,  we  have  computed 
the  mean  switching  time  as  the  noise  drives  the  system  out  of  the  basin  of  attraction  of  x a  to  x  <  xs-  In  each  of 
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Figure  5:  Log  of  escape  time  vs.  delay  (7  =  0.1).  The  bold  solid  curve  represents  the  theory  from  Eqs.  (30) 
and  (35),  and  it  is  shifted  vertically  by  -0.07  and  plotted  as  the  light  solid  curve  to  show  agreement  with 
the  data  in  rate  of  change.  The  data  points  are  the  mean  values  taken  over  2000  simulations. 

Figs.  4-6  the  solid  curve  represents  the  theory,  while  the  data  points  are  the  mean  values  taken  over  2000  simulations. 
To  insure  that  trajectory  doesn’t  drift  back  to  the  original  basin  we  define  the  switching  time  to  be  when  x  <  —0.01. 
To  find  the  theoretical  switching  time,  we  approximate  the  switching  rate  using  Eq.  (30),  with  the  minimum  7 Z  given 
by  Eq.  (35),  and  the  switching  time  is  Ts  =  1/W.  Here  we  use  the  shift  given  in  the  case  when  the  delay  is  small 
[40]. 

In  Fig.  4,  we  plot  the  switching  time  for  the  multiplicative  noise  case  as  a  function  of  the  inverse  of  the  noise 
intensity.  As  expected,  when  the  noise  intensity  is  increased  the  switching  time  decreases.  In  the  top  curve  the  delay 
is  set  to  t  =  0.1,  while  in  the  bottom  curve  the  delay  is  r  =  0.5.  As  before,  the  effect  of  the  delay  is  to  decrease  the 
switching  time  or,  in  other  words,  the  delay  increases  the  probability  of  switching. 

In  Fig.  5,  we  fix  the  noise  intensity  at  three  different  values  and  vary  the  delay  r.  There  is  a  small  shift  offsetting 
the  data  and  the  theory.  In  this  figure,  the  theory  was  shifted  vertically  by  -0.07  and  replotted  by  the  light  black 
solid  curve  to  show  agreement  with  the  data  in  rate  of  change.  Our  theory  based  on  small  r  maintains  good  fidelity 
up  almost  r  =  0.5;  for  t  >  0.5  the  theory  and  simulation  results  begin  to  diverge. 

In  Fig.  6,  the  delay  is  fixed  but  the  strength  of  the  delay  term  7  is  varied  from  zero  to  almost  one.  We  note 
that  the  theoretical  curve  deviates  from  the  simulation  results  as  7  gets  close  to  1;  this  is  due  to  a  singularity  in  the 
expression  for  the  shift  constant  c. 
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Figure  6:  Log  of  escape  time  vs.  dissipation  (r  =  0.1).  The  solid  curve  represents  the  theory  from  Eqs.  (30) 
and  (35),  while  the  data  points  are  the  mean  values  taken  over  2000  simulations. 

V.  CONCLUSIONS 

We  considered  the  problem  of  noise  induced  escape  from  a  basin  of  attraction  in  stochastic  systems  having 
differential  delay.  The  formulation  for  determining  the  escape  was  posed  as  a  variational  problem  along  an  optimal 
path,  which  maximized  the  probability  of  escape.  Since  noise  escape  from  a  basin  of  attraction  is  a  rare  event,  the 
optimized  solution  to  the  variational  problem  holds  in  the  tail  of  the  distribution.  The  scaling  results  derived  are 
therefore  based  on  using  the  optimized  path  to  escape. 

In  contrast  to  many  other  approaches  which  consider  Gaussian  processes,  the  approach  considered  here  is  quite 
general  in  that  it  allows  for  any  general  noise,  such  as  Poisson,  Schott,  or  colored  noise.  It  generalizes  the  problems 
of  non-Gaussian  noise  induced  escape  considered  in  [35,  38,  41]  by  extending  the  formulation  to  include  delay. 

We  applied  the  theory  to  scalar  problems  with  both  additive  and  multiplicative  noise  to  explicitly  compute  the 
solution  of  the  optimal  path  to  escape.  With  the  optimal  solution  in  hand,  we  derived  relations  describing  how  the 
escape  rate  scales  with  respect  to  various  parameters,  such  as  delay,  dissipation,  and  noise  intensity.  The  perturbation 
theory  used  in  the  analysis  was  based  upon  a  Melnikov  approach.  Here,  the  delay  itself  does  not  need  to  be  small, 
but  the  escape  path  with  delay  is  assumed  to  remain  close  to  the  path  with  no  delay.  Then,  a  first  order  correction 
to  the  integral  along  the  optimal  path  may  be  computed. 

The  escape  rate  results  found  analytically  were  compared  in  each  case  to  numerically  simulated  solutions.  Even 
though  a  small  delay  approximation  is  used  to  explicitly  compute  the  escape  exponents,  we  find  excellent  agreement 
between  theory  and  simulation  for  large  delay,  even  when  they  are  of  order  unity.  Specifically,  in  some  cases  for  the 
examples  considered  here,  the  escape  rate  agreement  between  theory  and  numerics  holds  for  delays  much  longer  than 
the  relaxation  time  of  the  system. 

Finally,  in  comparing  the  numerical  simulations  to  the  theory,  we  show  where  theory  and  numerical  simulations 
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diverge.  We  expect  that  the  approach  considered  here  will  shed  light  on  many  other  escape  processes  where  the 
stochastic  forcing  renders  the  problems  non-Markovian  in  nature,  as  well  as  many  other  applications  physical  and 
biological  experiments. 
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VI.  APPENDIX 


A.  Lagrangian  formulation 


The  Lagrangian,  L,  and  its  corresponding  second  order  equations  of  motion  for  a  stochastic  differential  delay 
equation  are  given  by: 


L(x,x,xT)  =  \[x{t)  -  F(x{t),x(t-r)]2, 

dL .  .  d  dL .  .  dL  , 

~  dtdA{x’X’XT)  =  ~d^{x-T’X-^x)- 


(36) 

(37) 


B.  Derivation  of  Delayed  Hamilton’s  Equations  of  Motion 

When  considering  the  Hamiltonian  with  delay,  new  equations  will  describe  the  momenta  due  to  the  addition  of 
a  time  delay.  This  is  now  derived,  and  shows  how  the  advanced  time  terms  arise. 

We  begin  by  assuming  the  Lagrangian,  L(x,x,xT),  of  the  systems  is  defined  as  a  function  of  state  vector  x, 
velocity  x ,  and  delayed  state  xT  =  x(t  —  r).  The  action  is  defined  as 


OO 

i in  J 


J[x,x\  =  min  /  L(x,x,xT)dt. 


(38) 


The  minimization  is  taken  over  all  paths  perturbed  from  the  true  minimum,  assuming  it  exists.  Defining  the 
Hamiltonian  using  the  usual  Legendre  transformation,  we  have  the  following  minimization  problem: 


OO 

JH[x,p]=min  J  (px  —  H(x,p,xT))dt, 


(39) 


where  p  =  x ,  x,xr).  We  assume  sufficient  hypotheses  so  that  the  inverse  function  theorem  applies.  We  let  the 
optimum  path  by  given  by  x  =  x0(t),p  =  pa(t),  and  consider  smooth  perturbations  from  the  path  as  T](t),^{t), 
respectively  so  that  limt^±ooTi{t)  =  0,  and  similarly  for  £(t).  That  is,  the  perturbed  paths  satisfy  the  natural 
boundary  conditions. 

Consider  the  variation  from  the  path: 


SJH[x0,Po\  =  J[xa  +  ri,p0  +  £}  -  J[xo,p0\ 

OO 

/dH  dH 

[pop  +  X o£  -  -T^(x0,po,x0r)p  -  —(x0,P0,x0t)£ 


Using  the  following  identities, 


-^-(x0,p0,x0T)pT]dt  +  O[\p\2  +  |£|2]. 


Popdt  =  -  popd-t, 


(40) 


(41) 


J  ^-{x0,p0,x0T)pTdt  =  [  ^-(x0(t  +  r),p0(t  +  T),x0(t))p(t)dt, 


dxT 


we  have  that 
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oo 

/dH 

[{Xo - —(x0,po,x0T))£(t) 


dH  dH 

-(po  +  (xo,po,  Xo-r)  +  Q^-(xo(t  +  T),po(t  +  T),xo(t)))ri(t)]dt 

+oi\v\2  +  m2]. 


Since  the  perturbations  to  the  path  are  arbitrary,  we  get  Hamilton’s  equations  of  motion  for  a  delayed  stochastic 
system  given  by  Eq.  (17). 


Xo 

Po 


—  (x0,P0,x0t) 

dH 

--q^-(xo,Po,XOt) 


dH 

dxT 


(x o(t  +  r),p0(t  +  r),  x0 (i)). 


(42) 


